<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/xhtml;charset=UTF-8"/>
<meta http-equiv="X-UA-Compatible" content="IE=9"/>
<meta name="generator" content="Doxygen 1.9.1"/>
<meta name="viewport" content="width=device-width, initial-scale=1"/>
<title>Eigen: Solving linear least squares systems</title>
<link href="tabs.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="jquery.js"></script>
<script type="text/javascript" src="dynsections.js"></script>
<link href="navtree.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="resize.js"></script>
<script type="text/javascript" src="navtreedata.js"></script>
<script type="text/javascript" src="navtree.js"></script>
<link href="search/search.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="search/searchdata.js"></script>
<script type="text/javascript" src="search/search.js"></script>
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
  $(document).ready(function() { init_search(); });
/* @license-end */
</script>
<script type="text/x-mathjax-config">
  MathJax.Hub.Config({
    extensions: ["tex2jax.js", "TeX/AMSmath.js", "TeX/AMSsymbols.js"],
    jax: ["input/TeX","output/HTML-CSS"],
});
</script>
<script type="text/javascript" async="async" src="https://cdn.mathjax.org/mathjax/latest/MathJax.js"></script>
<link href="doxygen.css" rel="stylesheet" type="text/css" />
<link href="eigendoxy.css" rel="stylesheet" type="text/css">
<!--  -->
<script type="text/javascript" src="eigen_navtree_hacks.js"></script>
</head>
<body>
<div id="top"><!-- do not remove this div, it is closed by doxygen! -->
<div id="titlearea">
<table cellspacing="0" cellpadding="0">
 <tbody>
 <tr style="height: 56px;">
  <td id="projectlogo"><img alt="Logo" src="Eigen_Silly_Professor_64x64.png"/></td>
  <td id="projectalign" style="padding-left: 0.5em;">
   <div id="projectname"><a href="http://eigen.tuxfamily.org">Eigen</a>
   &#160;<span id="projectnumber">3.4.90 (git rev 67eeba6e720c5745abc77ae6c92ce0a44aa7b7ae)</span>
   </div>
  </td>
   <td>        <div id="MSearchBox" class="MSearchBoxInactive">
        <span class="left">
          <img id="MSearchSelect" src="search/mag_sel.svg"
               onmouseover="return searchBox.OnSearchSelectShow()"
               onmouseout="return searchBox.OnSearchSelectHide()"
               alt=""/>
          <input type="text" id="MSearchField" value="Search" accesskey="S"
               onfocus="searchBox.OnSearchFieldFocus(true)" 
               onblur="searchBox.OnSearchFieldFocus(false)" 
               onkeyup="searchBox.OnSearchFieldChange(event)"/>
          </span><span class="right">
            <a id="MSearchClose" href="javascript:searchBox.CloseResultsWindow()"><img id="MSearchCloseImg" border="0" src="search/close.svg" alt=""/></a>
          </span>
        </div>
</td>
 </tr>
 </tbody>
</table>
</div>
<!-- end header part -->
<!-- Generated by Doxygen 1.9.1 -->
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
var searchBox = new SearchBox("searchBox", "search",false,'Search','.html');
/* @license-end */
</script>
</div><!-- top -->
<div id="side-nav" class="ui-resizable side-nav-resizable">
  <div id="nav-tree">
    <div id="nav-tree-contents">
      <div id="nav-sync" class="sync"></div>
    </div>
  </div>
  <div id="splitbar" style="-moz-user-select:none;" 
       class="ui-resizable-handle">
  </div>
</div>
<script type="text/javascript">
/* @license magnet:?xt=urn:btih:cf05388f2679ee054f2beb29a391d25f4e673ac3&amp;dn=gpl-2.0.txt GPL-v2 */
$(document).ready(function(){initNavTree('group__LeastSquares.html',''); initResizable(); });
/* @license-end */
</script>
<div id="doc-content">
<!-- window showing the filter options -->
<div id="MSearchSelectWindow"
     onmouseover="return searchBox.OnSearchSelectShow()"
     onmouseout="return searchBox.OnSearchSelectHide()"
     onkeydown="return searchBox.OnSearchSelectKey(event)">
</div>

<!-- iframe showing the search results (closed by default) -->
<div id="MSearchResultsWindow">
<iframe src="javascript:void(0)" frameborder="0" 
        name="MSearchResults" id="MSearchResults">
</iframe>
</div>

<div class="header">
  <div class="headertitle">
<div class="title">Solving linear least squares systems<div class="ingroups"><a class="el" href="group__DenseLinearSolvers__chapter.html">Dense linear problems and decompositions</a></div></div>  </div>
</div><!--header-->
<div class="contents">
<p>This page describes how to solve linear least squares systems using Eigen. An overdetermined system of equations, say <em>Ax</em> = <em>b</em>, has no solutions. In this case, it makes sense to search for the vector <em>x</em> which is closest to being a solution, in the sense that the difference <em>Ax</em> - <em>b</em> is as small as possible. This <em>x</em> is called the least square solution (if the Euclidean norm is used).</p>
<p>The three methods discussed on this page are the SVD decomposition, the QR decomposition and normal equations. Of these, the SVD decomposition is generally the most accurate but the slowest, normal equations is the fastest but least accurate, and the QR decomposition is in between.</p>
<h1><a class="anchor" id="LeastSquaresSVD"></a>
Using the SVD decomposition</h1>
<p>The <a class="el" href="classEigen_1_1SVDBase.html#ab28499936c0764fe5b56b9f4de701e26">solve() </a> method in the <a class="el" href="classEigen_1_1BDCSVD.html" title="class Bidiagonal Divide and Conquer SVD">BDCSVD</a> class can be directly used to solve linear squares systems. It is not enough to compute only the singular values (the default for this class); you also need the singular vectors but the thin SVD decomposition suffices for computing least squares solutions:</p>
<table class="example">
<tr>
<th>Example:</th><th>Output: </th></tr>
<tr>
<td><div class="fragment"><div class="line"><span class="preprocessor">#include &lt;iostream&gt;</span></div>
<div class="line"><span class="preprocessor">#include &lt;Eigen/Dense&gt;</span></div>
<div class="line"> </div>
<div class="line"><span class="keywordtype">int</span> main()</div>
<div class="line">{</div>
<div class="line">   <a class="code" href="classEigen_1_1Matrix.html">Eigen::MatrixXf</a> A = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">Eigen::MatrixXf::Random</a>(3, 2);</div>
<div class="line">   std::cout &lt;&lt; <span class="stringliteral">&quot;Here is the matrix A:\n&quot;</span> &lt;&lt; A &lt;&lt; std::endl;</div>
<div class="line">   <a class="code" href="classEigen_1_1Matrix.html">Eigen::VectorXf</a> b = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">Eigen::VectorXf::Random</a>(3);</div>
<div class="line">   std::cout &lt;&lt; <span class="stringliteral">&quot;Here is the right hand side b:\n&quot;</span> &lt;&lt; b &lt;&lt; std::endl;</div>
<div class="line">   std::cout &lt;&lt; <span class="stringliteral">&quot;The least-squares solution is:\n&quot;</span></div>
<div class="line">        &lt;&lt; A.template bdcSvd&lt;Eigen::ComputeThinU | Eigen::ComputeThinV&gt;().solve(b) &lt;&lt; std::endl;</div>
<div class="line">}</div>
<div class="ttc" id="aclassEigen_1_1DenseBase_html_ae814abb451b48ed872819192dc188c19"><div class="ttname"><a href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">Eigen::DenseBase::Random</a></div><div class="ttdeci">static const RandomReturnType Random()</div><div class="ttdef"><b>Definition:</b> Random.h:114</div></div>
<div class="ttc" id="aclassEigen_1_1Matrix_html"><div class="ttname"><a href="classEigen_1_1Matrix.html">Eigen::Matrix</a></div><div class="ttdoc">The matrix class, also used for vectors and row-vectors.</div><div class="ttdef"><b>Definition:</b> Matrix.h:182</div></div>
</div><!-- fragment -->  </td><td><pre class="fragment">Here is the matrix A:
  0.68  0.597
-0.211  0.823
 0.566 -0.605
Here is the right hand side b:
 -0.33
 0.536
-0.444
The least-squares solution is:
-0.67
0.314
</pre>   </td></tr>
</table>
<p>This is example from the page <a class="el" href="group__TutorialLinearAlgebra.html">Linear algebra and decompositions </a>. If you just need to solve the least squares problem, but are not interested in the SVD per se, a faster alternative method is <a class="el" href="classEigen_1_1CompleteOrthogonalDecomposition.html" title="Complete orthogonal decomposition (COD) of a matrix.">CompleteOrthogonalDecomposition</a>.</p>
<h1><a class="anchor" id="LeastSquaresQR"></a>
Using the QR decomposition</h1>
<p>The solve() method in QR decomposition classes also computes the least squares solution. There are three QR decomposition classes: <a class="el" href="classEigen_1_1HouseholderQR.html" title="Householder QR decomposition of a matrix.">HouseholderQR</a> (no pivoting, fast but unstable if your matrix is not rull rank), <a class="el" href="classEigen_1_1ColPivHouseholderQR.html" title="Householder rank-revealing QR decomposition of a matrix with column-pivoting.">ColPivHouseholderQR</a> (column pivoting, thus a bit slower but more stable) and <a class="el" href="classEigen_1_1FullPivHouseholderQR.html" title="Householder rank-revealing QR decomposition of a matrix with full pivoting.">FullPivHouseholderQR</a> (full pivoting, so slowest and slightly more stable than <a class="el" href="classEigen_1_1ColPivHouseholderQR.html" title="Householder rank-revealing QR decomposition of a matrix with column-pivoting.">ColPivHouseholderQR</a>). Here is an example with column pivoting:</p>
<table class="example">
<tr>
<th>Example:</th><th>Output: </th></tr>
<tr>
<td><div class="fragment"><div class="line"><a class="code" href="group__matrixtypedefs.html#ga731599f782380312960376c43450eb48">MatrixXf</a> A = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">MatrixXf::Random</a>(3, 2);</div>
<div class="line"><a class="code" href="group__matrixtypedefs.html#ga8028d921d43acd5605eabad41c254ef2">VectorXf</a> b = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">VectorXf::Random</a>(3);</div>
<div class="line">cout &lt;&lt; <span class="stringliteral">&quot;The solution using the QR decomposition is:\n&quot;</span></div>
<div class="line">     &lt;&lt; A.<a class="code" href="classEigen_1_1MatrixBase.html#adee8c19c833245bbb00a591dce68e8a4">colPivHouseholderQr</a>().solve(b) &lt;&lt; endl;</div>
<div class="ttc" id="aclassEigen_1_1MatrixBase_html_adee8c19c833245bbb00a591dce68e8a4"><div class="ttname"><a href="classEigen_1_1MatrixBase.html#adee8c19c833245bbb00a591dce68e8a4">Eigen::MatrixBase::colPivHouseholderQr</a></div><div class="ttdeci">const ColPivHouseholderQR&lt; PlainObject &gt; colPivHouseholderQr() const</div><div class="ttdef"><b>Definition:</b> ColPivHouseholderQR.h:664</div></div>
<div class="ttc" id="agroup__matrixtypedefs_html_ga731599f782380312960376c43450eb48"><div class="ttname"><a href="group__matrixtypedefs.html#ga731599f782380312960376c43450eb48">Eigen::MatrixXf</a></div><div class="ttdeci">Matrix&lt; float, Dynamic, Dynamic &gt; MatrixXf</div><div class="ttdoc">Dynamic×Dynamic matrix of type float.</div><div class="ttdef"><b>Definition:</b> Matrix.h:500</div></div>
<div class="ttc" id="agroup__matrixtypedefs_html_ga8028d921d43acd5605eabad41c254ef2"><div class="ttname"><a href="group__matrixtypedefs.html#ga8028d921d43acd5605eabad41c254ef2">Eigen::VectorXf</a></div><div class="ttdeci">Matrix&lt; float, Dynamic, 1 &gt; VectorXf</div><div class="ttdoc">Dynamic×1 vector of type float.</div><div class="ttdef"><b>Definition:</b> Matrix.h:500</div></div>
</div><!-- fragment -->  </td><td><pre class="fragment">The solution using the QR decomposition is:
-0.67
0.314
</pre>   </td></tr>
</table>
<h1><a class="anchor" id="LeastSquaresNormalEquations"></a>
Using normal equations</h1>
<p>Finding the least squares solution of <em>Ax</em> = <em>b</em> is equivalent to solving the normal equation <em>A</em><sup>T</sup><em>Ax</em> = <em>A</em><sup>T</sup><em>b</em>. This leads to the following code</p>
<table class="example">
<tr>
<th>Example:</th><th>Output: </th></tr>
<tr>
<td><div class="fragment"><div class="line"><a class="code" href="group__matrixtypedefs.html#ga731599f782380312960376c43450eb48">MatrixXf</a> A = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">MatrixXf::Random</a>(3, 2);</div>
<div class="line"><a class="code" href="group__matrixtypedefs.html#ga8028d921d43acd5605eabad41c254ef2">VectorXf</a> b = <a class="code" href="classEigen_1_1DenseBase.html#ae814abb451b48ed872819192dc188c19">VectorXf::Random</a>(3);</div>
<div class="line">cout &lt;&lt; <span class="stringliteral">&quot;The solution using normal equations is:\n&quot;</span></div>
<div class="line">     &lt;&lt; (A.transpose() * A).ldlt().solve(A.transpose() * b) &lt;&lt; endl;</div>
</div><!-- fragment -->  </td><td><pre class="fragment">The solution using normal equations is:
-0.67
0.314
</pre>   </td></tr>
</table>
<p>This method is usually the fastest, especially when <em>A</em> is "tall and skinny". However, if the matrix <em>A</em> is even mildly ill-conditioned, this is not a good method, because the condition number of <em>A</em><sup>T</sup><em>A</em> is the square of the condition number of <em>A</em>. This means that you lose roughly twice as many digits of accuracy using the normal equation, compared to the more stable methods mentioned above. </p>
</div><!-- contents -->
</div><!-- doc-content -->
<!-- start footer part -->
<div id="nav-path" class="navpath"><!-- id is needed for treeview function! -->
  <ul>
    <li class="footer">Generated on Thu Apr 21 2022 13:07:55 for Eigen by
    <a href="http://www.doxygen.org/index.html">
    <img class="footer" src="doxygen.png" alt="doxygen"/></a> 1.9.1 </li>
  </ul>
</div>
</body>
</html>
